Electronic-enthalpy functional for finite systems under pressure 
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We introduce the notion of electronic enthalpy for first-principles structural and dynamical cal- 
culations of finite systems under pressure. An external pressure field is allowed to act directly on 
the electronic structure of the system studied via the ground-state minimization of the functional 
E + PV q , where V q is the quantum volume enclosed by a charge isosurface. The Hellmann-Feynman 
theorem applies, and assures that the ionic equations of motion follow an isoenthalpic dynamics. No 
pressurizing medium is explicitly required, while coatings of environmental ions or ligands can be 
introduced if chemically relevant. We apply this novel approach to the study of group-IV nanopar- 
ticles during a shock wave, highlighting the significant differences in the plastic or elastic response 
of the diamond cage under load, and their potential use as novel nanostructured impact-absorbing 
materials. 
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The study of nanoparticles under pressure is rapidly 
acquiring great scientific and technological interest, since 
it allows to explore structural arrangements and phase 
transitions under different kinetic and thermodynamic 
conditions than those of bulk solids. While experiments 
on the compression of nanoparticles have been performed 
for more than a decade 0J3 l _P, 0, Q, and encompass 
broad classes of materials [fj, [7|, computational studies 
have been hampered by the conceptual difficulty arising 
in applying classical or electronic-structure methods to 
the study of finite systems under pressure. Extended 
systems can be studied using variable-cell isoenthalpic 
dynamics H |j| 0, |n|; however, these methods do not 
carry over to the case of finite systems, unless the en- 
vironment (i.e. the pressuri zing m edium) is introduced 
explicitly in the calculations |12 . ll3lll3 |. A full quantum- 
mechanical treatment becomes then too costly in all but 
the smallest cases, due to the large number of atoms or 
molecules that are needed to reproduce the pressuriz- 
ing environment; in addition, equilibration of the system 
at a given pressure requires extensive simulations aver- 
aging over many collision events. To solve the first of 
these problems, a mixed quantum-mechanics/molecular- 
mechanics approach has been proposed [I2I where 
the quantum fragment is compressed by a pressurizing 
medium described by classical force fields. 

In this Letter we present a novel and general approach 
to this problem, and show how pressure can be applied 
to a finite system without the need to introduce a pres- 
surizing medium. We introduce here, as thermodynamic 
functional describing a finite system under pressure, the 
electronic enthalpy 



chosen, and can be straightforwardly computed as 



H = E + PV a 



(1) 



where E is the internal energy of the system, P is the 
desired pressure and V q is the "quantum volume" occu- 
pied by the electronic charge density. This quantity is 
well defined once the threshold for a density isosurface is 



V q = I dr #(p(r) 



a) 



(2) 



where p(r) is the electronic density and § is a step func- 
tion at the threshold value a. For computational conve- 
nience we smooth the § to defined as the integral of a 
normalized Gaussian of width a = a/3 |l5|| . 

In density-functional theory, the total potential acting 
on the electrons will be given by the functional derivative 
of the electronic enthalpy in with respect to the charge 
density p(r). This results in an additional contribution 
deriving from the term PV q which is simply 



-0(r)-a) 2 /2<r 2 



(3) 



In the self-consistent solution of the Kohn-Sham equa- 
tions the $>v potential drives the evolution of the elec- 
tronic ground-state to the minimum of the electronic en- 
thalpy The Hellmann-Feynman theorem applies, 
and thus the ionic forces take implicitly into account 
the contribution of the external load, transferring the 
effects of the compression directly to the ionic relaxation 
and dynamics. A Lagrangian formulation (e.g. for first- 
principles Car-Parrinello molecular dynamics 

H3) 

also 

follows directly. 

The introduction of an electronic-enthalpy functional 
has several conceptual and practical advantages: 1) The 
pressure field acts directly on the electrons; thus, the 
compressibility of the system is properly dominated by 
electrostatic and Pauli-principle effects. Also, if pressure 
were transferred directly to the ionic nuclei by an ex- 
ternal classical force field, electrons could ultimately be 
squeezed out of the system. 2) Isoenthalpic relaxations 
and dynamics do not require extensive equilibration with 
a pressurizing medium; a medium can always be intro- 
duced as an environmental coating, if chemically rele- 
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FIG. 1: The dependence of the innermost interlayer dis- 
tance on pressure for a 16-layer silicon slab relaxed using 
the electronic-enthalpy functional (solid curve); results from 
the Murnaghan equation of state for a 16-layer bulk super- 
cell (dotted) are shown for comparison. The inset shows the 
results obtained for a 8-layer slab (solid), a 8-layer supercell 
(dotted), and a 8-layer slab compressed via the explicit intro- 
duction of argon as a pressurizing medium (dashed). 



vant. 3) In the thermodynamic limit, the classical macro- 
scopic volume and enthalpy are recovered. 4) Implemen- 
tation in electronic-structure codes is straightforward, 
and can be easily extended to other wavefunction-based 
approaches such as Hartree-Fock or variational quantum 
Monte Carlo. (The results presented here have all been 
obtained using first-principles Car-Parrinello molecular 
dynamics as implemented in the public domain 
^-ESPRESSO package [13.) 

Other developments and applications stemming from 
this formulation can be envisioned. A surface tension 
term H — > H + jS can be introduced, where the surface 
area is obtained by finite differences via 



S 



(4) 

and A determines the width of the discretization. This 
latter contribution to the energy functional could be used 
to characterize surface tension in metallic nanoparticles, 
confinement effects in electron bubbles (e.g., for electrons 
in superfluid helium 0) or cavitation effects in solva- 
tion models [l^. Finally, we stress that although classi- 
cal definitions of volume can be introduced in atomistic 
simulations [2(1 |2l|, a natural measure based on the elec- 
tron density can only emerge in a quantum-mechanical 
formalism. A validation of the present approach can be 
provided not only in the thermodynamical limit (where 
the quantum volume V q becomes indistinguishable from 
the classical limit), but also for finite, small systems. For 
this purpose we studied the effect of pressure on thin sil- 
icon slabs. This reference system allows for a direct com- 
parison with bulk results without introducing additional 



surface tension effects. We plot in Fig. ^the average dis- 
tance between the innermost layers in two Si(OOl) slabs 
of 8 and 16 layers, as a function of the pressure parameter 
P entering our functional (Q. A remarkable agreement 
is found with the results obtained when calculating the 
physical pressure [^ in a bulk silicon system having the 
same geometry and compressed along the [001] direction. 
For the 8-layer slab a full quantum-mechanical calcula- 
tion taking into account explicitly argon as a pressurizing 
medium has been performed, delivering an even closer 
agreement (at least up to 20 GPa) to our electronic- 
enthalpy results. This agreement clearly demonstrates 
that the minimization of the quantum enthalpy correctly 
accounts for the effects of pressure; in addition, agree- 
ment is reached rapidly with respect to system size |23|. 

Our case studies are group- IV nanoparticles M35H36 
during a shock wave (M corresponds to Si, Ge and C). 
The initial configurations, carved from bulk diamond, 
have all dangling bonds saturated with hydrogens. The 
structures were relaxed and thermalized for 2-3 ps at 300 
K and zero pressure, at which point the ionic thermo- 
stat was switched off (no electronic thermostat is needed 
in these simulations). We modeled the adiabatic com- 
pression characterizing a shock wave with a linear in- 
crease in hydrostatic pressure from to 40 GPa for 1 
ps, followed by 0.75 ps at constant pressure and by a 
linear decrease to GPa in 1 ps. The system was fi- 
nally allowed to evolve at zero pressure for another 8- 
10 ps. The structural transformations occurring dur- 
ing the shock are highlighted in Figs. [5] and 03 where 
the time evolution of the average bond length, quan- 
tum volume, and average coordination of the ions are 
shown for the three nanoparticles. The response of the 
Si nanoparticle (see Fig. |3J) is remarkable: a transition is 
induced ^0.6 ps into the compression (P ~ 26-28 GPa), 
at which point the system amorphizes. In this process 
the second shell of neighbors overlaps and merges with 
the first, increasing the coordination of the ions (see Fig. 
|3J) from 3.0 - due to the undercoordinated surface atoms 
- to ^5.5. The compressed bonds recover, albeit only 
on average, the zero-pressure equilibrium length; volume 
contraction is achieved by a collapse of the coordination 
shells rather than a shortening of the nearest bonds. The 
resulting amorphous configuration remains stable after 
pressure is released, and no recrystallization is detected 
in the following 8-10 ps. A structural relaxation con- 
firms that the amorphous nanoparticle corresponds to a 
well-defined local minimum for the potential energy sur- 
face; this minimum is 11.15 eV (0.32 eV per silicon atom) 
higher than the corresponding initial tetrahedral config- 
uration. These results confirm the existence of an amor- 
phous phase for silicon nanoparticles, as reported in Ref. 
P and also found in the simulations of Refs. 0, Il3j |. 
Germanium and carbon at these sizes display instead an 
elastic response (see Fig. |2J) and the original tetrahedral 
structure is recovered at the end of the simulation. Due 
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FIG. 2: Time evolution of the average bond length (solid 
line) and of the quantum volume (dashed line) for the three 
nanoparticles during the shock wave. 



to its very stiff bonds, the carbon nanoparticle maintains 
its tetrahedral coordination even at peak pressure, while 
germanium is noticeably softer. In both cases, however, 
the average coordination of the ions is barely affected (see 
Fig. |3J) and the first two peaks in the pair-correlation 
function remain distinct at all times during the shock. 
The differences observed in the structural behavior are 
also reflected in the thermodynamics of the compression 
process. In Fig. 0] we show the mechanical work ex- 
erted on the nanoparticles by the compression wave; in 
our approach this can be straightforwardly determined 
by j PdV, whereas the fluctuations of an explicit pres- 
surizing medium would hinder the volume changes of the 
quantum fragment. It is immediately evident that the 
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FIG. 3: Time evolution of the average ionic coordination for Si 
(solid line), Ge (dashed line) and C (dot-dashed line) nanopar- 
ticles during the shock wave. 



Si nanoparticle absorbs a great amount of mechanical 
energy: a portion of this energy is converted into ki- 
netic energy of the ions, thanks to enhanced anharmonic 
coupling between acoustic and optical phonons, while a 
larger part goes into the plastic deformation and tran- 
sition into the metastable amorphous phase. Both the 
diamond nanoparticle and (to a lesser extent) the ger- 
manium one release instead elastically all the mechanical 
energy absorbed during the shock and stored in the com- 
pression of the covalent bonds. 

Kinetic trapping for small Si nanoparticles has been 
demonstrated in Ref. Q, with an amorphous phase ap- 
pearing during pressure release below 5 GPa. In agree- 
ment with Ref. 0, we find an enhancement with re- 
spect to the bulk for the pressure transition during the 
upstroke phase; our nanoparticle falls directly in the 
amorphous phase, without an intermediate Imma or j3- 
Sn phase. Pressure-induced amorphization was observed 
when studying porous silicon films containing diamond- 
structure nanocrystals under load |24) : amorphization 
was also found in bulk first-principles simulations [jjlEol 
in the 11-15 GPa range. The importance of the rate of 
load release on the structural transformation in Si and 
Ge was highlighted in Ref. 27]. We have verified that 
changing the speed of the shock wave by a factor of 1/15, 
1/2, or the initial temperature of the system from 300 K 
to 100 K, did not change the amorphization pathway or 
the transition pressure. We also carefully confirmed that 
the nanoparticle never becomes liquid. 

Finally, we examined the effects that the threshold 
parameter a defining the quantum volume has on the 
isoenthalpic dynamics. The silicon simulation was re- 
peated using a = 0.0005 electrons/ ( a. u.) 3 (a = 0.0002 
electrons/ (a. u.) 3 was used for all previous simulations). 
While at GPa this different choice leads to a volume 
difference of ~ 10%, in the course of the compression the 
two definitions rapidly converge to the same value; the 
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FIG. 4: Work J PdV done by the external pressure on the 
nanoparticles during the shock wave (solid line) , and instanta- 
neous kinetic energy of the ions (dashed line, zero corresponds 
to the 300 K average). 

soft "electronic skin" involved when using a lower thresh- 
old contracts quickly, and the system then opposes the 
same resistance to the external load as in the case of a 
higher threshold. More importantly, the energy absorbed 
is shown to be the same independently of our choice of 
threshold; the outer "electronic skin" releases elastically 
at the end of the run any energy stored at the beginning 
of the simulation. 

In conclusion, we have introduced a novel approach 
in the study of finite systems under pressure, based 
on the quantum-mechanical definition of an electronic- 
enthalpy functional. This method is straightforward to 
implement in static and dynamical first-principles cal- 
culations, and opens the study of the phase stability 
and kinetics of finite quantum systems to extensive first- 
principles simulations. Its accuracy and broad indepen- 
dence on the isosurface threshold have also been estab- 
lished. We have studied the behavior of semiconducting 
nanoparticles (M35H36, with M = Si, Ge, and C) dur- 
ing a hydrostatic shock wave, highlighting the different 
response of the three materials. Small silicon nanopar- 
ticles undergo plastic amorphization at pressures larger 
than those required to amorphize a bulk system, while 
carbon nanoparticles (and to a lesser extent germanium 
ones) elastically recover their original tetrahedral struc- 
ture [28j. A large amount of mechanical energy is ab- 
sorbed by silicon nanoparticles as internal energy of their 
amorphous phase and kinetic energy of the ions. This 
energy will be released in longer timescales than those 
of a traveling shock wave, hinting at the possibility of 
designing a nanostructured impact-absorbing material. 
It is also envisioned that the transition pressures, decay 
timescales, and the overall mechanical response could be 
tuned by an appropriate choice of materials, composi- 
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